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Abstract 

An initial stress within a solid can arise to support external loads or from pro¬ 
cesses such as thermal expansion in inert matter or growth and remodelling in living 
materials. For this reason it is useful to develop a mechanical framework of initially 
stressed solids irrespective of how this stress formed. An ideal way to do this is to 
write the free energy density dt = in terms of initial stress r and the elastic 

deformation gradient F. In this paper we present a new constitutive condition for 
initially stressed materials, which we call the initial stress symmetry (ISS). We focus 
on two consequences of this symmetry. First we examine how ISS restricts the free 
energy density dt = d'(F,T) and present two examples of d'(F,T) that satisfy ISS. 
Second we show that the initial stress can be derived from the Cauchy stress and the 
elastic deformation gradient. To illustrate we take an example from biomechanics and 
calculate the optimal Cauchy stress within an artery subjected to internal pressure. 

We then use ISS to derive the optimal target residual stress for the material to achieve 
after remodelling. 

Keywords: residual stress, initial stress, biomechanics, elasticity, constitutive equations 


1 Introduction 


When all loads are removed a body can still hold a signihcant amount of internal stress, 
called the residual stress. In manufacturing, residual stress has long been noted to be 
detrimental to, or enhance, the performance of a material. For biological tissues, residual 
stress is used to self-regulate stress and strain, and ultimately preserve ideal mechanical 
conditions for the tissue (Fung, 1991 Holzapfel and Ogden, 2003). In geophysics, due to 
gravity the Earth has developed high initial stress within, which greatly influence the 
propagation of elastic waves. 

Here we use the term initial stress to broadly mean the internal stress of some reference 
configuration, irrespective of how the stress was formed or the boundary conditions. In 
this sense residual stress is a form of initial stress. 

The initial stress felt by any region of a material is due to the push and pull of the 
surrounding regions. If any region were to be cut out from the material, the stress on its 
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newly formed boundary would be zero, thus reducing the potential energy in the bulk. 
Based on this concept Roger developed constitutive laws for residually stressed materials, 
see for example Hogei (1986); Hogei (1993); Johnson and Hogei (1995) Roger showed 
that by taking this idea to its limit, and cutting the material into possibly an infinite 
number of disconnected regions, the material may be relieved of all of its internal stress in 
a configuration called the virtual stress-free state. From that configuration, a hyperelastic 
energy can be defined as a function of the strain from the virtual state to the current 
configuration. 

Though the use of this virtual state is technically sound, it leads to challenging cal¬ 
culations even for simple deformations and rarely yields analytic results, unless great 
simplifications are assumed. Moreover, the experimentally identification of the virtual 
state requires cutting the material, which is not always suitable, especially for living or¬ 
ganisms. Rowever, using a virtual stress-free reference is routinely seen as the only viable 
alternative, quoting Chuong and Fun ^ (1986) “To characterize the arterial wall or any 
other biological soft tissue, we need a stress-free state”. Conversely, We believe that by 
developing tools to work directly with initially stressed reference configurations, without 
the need of a stress-free state, will be very useful, specially in biomechanics. 

An ideal way to account for the initial stress would be to have a free energy density 
function T = T {F, r) written explicitly in terms of the deformation gradient F and the 
initial stress r, without any a priori restrictions. For the development of constitutive 
laws there is no need to distinguish between residual stresses and initial stresses, a view 
which is shared with Merodio et al (2013) The initial stress r could then be determined 


from elastic wave speeds (Shams and Ogden, 2014 Destrade and Ogden, 2012 Man and 
Lu, 1987) or by solving the linear equations of momentum balance. 

Shams et al. (2011)] and Shams (2010) worked towards a general framework for ini¬ 
tially stressed solids, while others have investigated the mechanics for some examples of 
^(F,t) (Merodio et ah, 2013; Merodio and Ogden, 2015). For a more complex geometry, 
Wang et al ] (2014')] found that including one residual stress invariant in the free energy den¬ 
sity was a simple way to model the effects of residual stress on the myocardium. Rowever, 
in general a major obstacle still remains: how to write 'k(F', r) in terms of the combined 
invariants of F and r? If the only source of anisotropy is due to the residual stress, then 
the free energy still depends on ten independent invariants. 

Johnson and Rogei |~(l995) developed representations for the Cauchy stress response 
a = q{F,T) in terms of F and r by assuming a stress free virtual state and numeri¬ 
cally inverting the residual stress-strain equation. Rowever, this approach often requires 
solving numerical nonlinear implicit equations. Johnson and Rogei (1995); Johnson and 
Rogei (1998) exemplified this approach for a material with virtual state composed by a 
Mooney-Rivlin strain energy. 

In this work, we introduce a new constitutive requirement on 'I'(i^, r) called the Initial 
Stress Symmetry (ISS). ISS restricts the constitutive form of the stress cr = ^{F,t) by 
providing a constitutive equation for the initial stress r = ^{F~^,<t) that must hold for 
every F and r. To our knowledge this symmetry has never been discussed before. 

The ISS also helps answer an important question: how much can the Cauchy stress be 
altered by adjusting the initial stress? From a modelling perspective, initial stress has been 
used to make the material more or less compliant (Johnson and Roger, 1998), to control the 
Poynting effect (Merodio et ah, 2013) or to maintain an ideal internal stress (Fung, 1991). 
Given a 'I'(i^, r) that satisfies ISS, then r = ^{F~^, a) suggests that for any choice of <t, 
there will exist an initial stress r that supports cr. 

The basic equations for an elastic material subject to initial stress are summarized 
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in Section and then the initial stress symmetry is first presented in Section 2.1 ISS 


is satisfied automatically if the initial stress is due to an elastic deformation of a stress 
free configuration, which we demonstrate in Appendix For this reason we develop in 
Section|^an example for that satisfies ISS, by deforming an incompressbile neo- 

Hookean material from a stress free virtual state. 

In Sectionj^we express ISS, in all generality, as nine scalar equations for an incompress¬ 
ible material, written in terms of two undetermined scalars p and Pr] and the invariants 
of F and r. This form of ISS makes it easier to select representations for ^{F,t). For 
example, in Section 4.1 we show, with a minor adjustament to the equations, how to use 
the scalar equations of ISS to deduce an example of 4'(F', r) for a compressibly material. 

It is commonly thought that arteries attempt to maintain a homogeneous stress gra¬ 
dient within their walls (Taber and Eggers, 1996). In Section]^ we calculate this optimal 
Cauchy stress for a simplified arterial wall, and then show how by using ISS we can calcu¬ 
late the residual stress that exactly supports this optimal Cauchy stress in SectionWe 
finally compare these results against the commonly used opening-angle method (Chuong 
and Fung, pT86 ). 


2 Initially stressed elastic materials 


A common approach to model the effect of residual stress is to consider a virtual stress free 
configuration B, from which the material is deformed and “glued” together to produce a 


residually stressed equilibrium state B. See Figure 
An elastic stored energy density T can then be c 
gradient F from B to the current configuration B so that 'k = 'k(F'). 


^for a diagram of all the configurations, 
efined as a function of the deformation 


O 

F 



Figure I: B is the current configuration with internal stress cr, while .6 is a reference 
configuration with internal stress r. The virtual stress-free state is a collection of 
configurations where the body is stress-free. 

In this paper we want to write the free energy density T as a function of the initial 
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stress T and of the deformation gradient F : B ^ B, so that 'I' = r), and B will not 

necessarily be an unloaded configuration. We feel that it is natural to consider that initial 
stress contributes to the potential energy stored by a material. One extreme example is 
the wapa tree, which has been know to burst open once cut, possibly causing injury, due 


to its immense level of residual stress (Detienne and Thiel, 1988). 

We assume that F is a purely elastic deformation, but we do not make any assumptions 
about the origins of the initial stress r, except that r affects the stored energy density 'k. 
Assuming that the body is incompressible, i.e. J = det F" = 1 at all times, the Cauchy 


stress tensor a reads (Guillou and Ogden, 2006): 


(O'k 


( 1 ) 


where p is the Lagrange multiplier associated with the constraint of incompressibility, I 
is the identity matrix and r is the initial stress. If the material is compressible then p is 
replaced by —2I^d^/dI^. Note we have and will omit the possible dependence of 'k on 

O O 

the position X € B for the sake of simplicity. For the body in the configuration B we have 
that F = I and cr = t; thus we require that 


r = ^{I,T)-pI, (2) 

where p is the value of p when F = I. We call the above equation the residual stress 
compatibility. 

The presence of residual stress generally leads to an anisotropy response of the material 

O 

in reference to B. Here we assume no other source of intrinsic anisotropic so that 'k can 
be written as a function of all the independent invariants generated by r and C = F'^F, 
the right Cauchy-Green deformation tensor. Following Shams et al.(2011), we take the 
following complete set of 10 independent invariants 


= tr C, h = \[{ll - tT{C% h = det C, 

= trr, Ira = = det T, 

Ji = tr(TC), J2 = tr(rC^), J3 = tr(T^C'), J4 = tr(T^C^). 


(3) 

(4) 

(5) 


The Cauchy stress Eq. Q can then be written as (Shams et al., [2011 ) 


a = 2^i^B + 2^/2(IiF - B^)-pI + 2 ^j,FtF^ + 2^j^F{tC + Ct)F^ 

+ 2^j^Ft^F^ + 2'I>j^F{t^C + Ct^)F^ ( 6 ) 


where B = FF^, and 'kjj,'kjj, are the partial derivatives of with 

respect to Ii, I2, Ji, J2, J3, Ji respectively. There are no partial derivative of ^ with respect 
to and Irg appearing in Q because r does not depend on F. 

The residual stress compatibility Eq. Q becomes 


2—— -|- 4—— 
9/1 dh 


-p = 0, 


dJi 


+ 4 


dJ2 


= 1 , 


94' 94' 

-h 2- 

9J3 9J4 


(7) 


Another important physical restriction that can be imposed is the strong-ellipticity 
condition, which is satisfied when the fourth-order tensor 


A 


Opiqj 


= J ^FpaFqjS 


92^ 

OFiadFjg’ 
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satisfies 


•A-opiqjnpnqmimj > 0 for every n,m G 


( 8 ) 


for compressible materials, while for incompressible materials the above need only hold for 
n ■ m = 0. Imposing SE implies that plane waves may propagate in every direction with 
a real valued speed (Truesdell, [1966 ), and other physically expected behaviour (Walton 


see Shams et al (2011) 


and Wilber, 2003). For a representation of Aopiqj in terms of the invariants of r and F 


The issue we address now is how to write 'h explicitly in terms of the invariants 0.0 
and (§? We advocate three criteria. The free energy density T should satify the initial 
stress compatibility 0 . Second, it should satisfy strong-ellipticity ( [s]) for a ll deformations 
in which the material is expected to be stable (Merodio and Ogden, 2003), and the third 
criterion we call the initial stress symmetry (ISS). 


2.1 Initial Stress Symmetry 

For convenience let the response function $ be denoted by 

d'it 

$(Fi,<T2,Pi) := i^i^(i^i,cr2) - Ipi, 


(9) 


for every Fi and cr 2 , where the argument on the right hand side is evaluated by taking 
the partial derivative of ^{F, ■) with respect to F. The scalar pi is undetermined if the 
material is incompressible and pi = — 2 / 3 T /3 with F replaced with i^i if the material is 
compressible. 

The ISS states that $ has no preferred reference configuration. Refering to Figure 

O 

if we take B as the reference configuration, then the Cauchy stress becomes 


cr = q{F,T,p). 


( 10 ) 


However, we can also take B as the reference configuration and B as the current configu¬ 
ration and therefore express the initial stress as 


T = $(F ,(T,Pr), 


( 11 ) 


for some scalar p,-. For a compressible material Pr = —213 ^with F replace with F 
In a more precise form, ISS can be stated as 


-1 


cr = q{F,T,p) and t = q{F Scr,p^), 


( 12 ) 


for every F and r, such that r = and det F = 1 for incompressible materials, andp and 
Pr are respectively given by p = p(cr, F, r) and Pr = p(r, F~^, a) for some scalar function 
p. The boundary conditions can determine p through its dependence on <t, and analogously 
for Pr- The ISS agrees with the initial stress compatibility, i.e. the condition Q, when 
we have ^{I,T,p) = r for every r. 

Another way to view the ISS symmetry is to assume that r and cr are due to the 
elastic deformation of a virtual stress-free state, which we demonstrate in Appendix 
Hoger emphasized many times that using the virtual stress-free state does not restrict 
how the residual stress was formed. The same can be said about ISS; only F is an elastic 
deformation. Moreover, the ISS is not restricted to elasticity and should hold for other 
constitutive equations such as those encountered in viscoelasticity and plasticity. 

One practi cal outcome from the ISS is to restrict the possible constitutive choices for 

In Section 4.1 we show that 'k = — 3) -|- |(Ji — 1^) proposed by Merodio et 
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al (2013) does not satisfy ISS, and we also use ISS to deduce an expression for 


for compressible materials. 

A second practical feature arising from ISS is that we can write the residual stress as 
a function of the Cauchy stress This proves useful when cr is known a priori, as it 

will be discussed in Section |5.2| where we determine cr from a homeostasis principal and 
then use Eq. © to derive r. The alternative of choosing the residual stress r first can 
be far more complicated. 

Before further developing the implications of ISS in Section]^ we will introduce below 
an example of stored energy ^{F, r) that satisfies ISS, stress compatibility Q and strong- 
ellipticity (Isj). 


3 Initially stressed neo-Hookean material 

Here we derive a simple constitutive equation for an initially stressed body by assuming 
that both T and cr arise from deforming an incompressible neo-Hookean material. The 


result will be an explicit representation of T in terms of F and r, given by Eq. (23) below, 
that automatically satisfies ISS, the stress compatibility ([^ and strong-ellipticity® . Both 


ISS and stress compatibility hold because r arises due to an elastic deformation from a 
stress-free state, see Appendix [A] and strong-ellipticity is satisfied because the material is 
a neo-Hookean solid ((^den, |1997[ ). 

Refering to Figure [ij as the material is stress-free in H, we have 


T = |(trC-3), 


(13) 


where /U > 0 is the constant shear modulus and C = F F. The Cauchy stress ([^ then 
becomes 

a = fiFp'^ -pi. (14) 

We can rewrite tr C by substituting F = FF and using the properties of the trace 

rp ^ O O O 

trC = tr(F F) = ^{F'^F'^FF) = tr(FF^F^F) = tr(SC), 
where B = FF^ and C = F^F, which leads to 


T = ^ 
2 


tr{BC) - 3 


(15) 


We can write B in terms of r by evaluating Eq. (14) at = I, cr = t, p = p and 
rearranging the term as follows 


^ O 

pB = T + pi. 

o o o 

To write p as a function of r we require that F be isochoric, so that det{pB) 
pi) = p?, which results in a cubic equation for p 

o3 o2 o o 

P +p lr^+ plr2 + It3 - P =0. 


(16) 
det(T -|- 

(17) 
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The real roots for p are given by, 


p = < 


C 1 T 3 + Ci— — Iri 

*rp T 

C 1 T 3 + Cl — — In 


T 2 < 


-T 


3/2 


1 < T 2 , 


Tn3/2 ^ rr ^ rp3/2 

-^1 S -^2 S -^1 ? 


(18) 


while if Ti = 0 then p = p — where 


Ti =1^^ - 3L 


27, 


^2) 


^2 — Iti + 2 ^^'^3 ^ 


T 3 =\ Jt^-t^-T2 , 


1 , Vs. 

Cl = "2 


(19) 

( 20 ) 


In terms of the eigenvalues ri, T 2 and of r, we can write 
27 1 

T 2 = + 2 [Vl “ "^ 2 ) - (r3 - Ti)] [{t2 - n) - (n - rs)] [(r3 - n) - (t2 - Tg)] , 

T’l = ^(n - rs)^ + ^("^1 - ”^ 3 )^ + ^V 2 - nf, (21) 

so that when Ti = 0 we have that ri = T 2 = T 3 . From the above we see that Ti > 0 and 
ji3/2 ^ ^ 21 12\ for anjf] Ti, r 2 , Ts G M. So if T 2 < 0 then clearly T 2 < else if 


T 2 > 0 then T 2 < IT 2 + /r^27/2| < Therefore T 2 < and so the condition for 


3/2 


the first case for p in Eq. (18) is always satisfied. We discard the second and third case 
for p because we expect T, and therefore p, to be continuous for every r G When 


Ti = T 2 = T 3 the only viable solution for j? is (18)i. If r moves into a region where (18)2 


or ( 18)3 becomes real, then p can not change from (18)i to ( 18)2 or ( 18)3 because it can 


be shown that (18)i does not equal (18)2 or ( 18)3 for any r G 


To represent tr(SC') in terms of the invariants of r and C, we multiply each side of 


Eq. (16) on the right with C and take the trace to get 


pti{BC) = tr(rC') +ptrC', 
which we use to write the free-energy density Eq.(|15[) 


as 


T = - 

2 


^ {ph + Ji - 3/x^ 


( 22 ) 


(23) 


with p given by Eq. (18)i. Note that in the absence of residual stress t = 0, p = p 


by ([TtI), Ji = h, and then T reduces to the classical neo-Hookean model, as expected. 
Equation (23) represents the general extension of the neo-Hookean strain energy function 


to a residually stressed material, resulting in a function of only five of the nine independent 
invariants of C and r. 

Substituting (23) in the Cauchy stress Eq. 0 we arrive at the constitutive relation 


<T = pB + FtF^ — pi. 


(24) 


‘http://math.stackexchange.com/questions/1255883/prove-that-a2b2c2-geq-2a-b2b-c2a-c21- 

3/1255907#1255907 
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Taking B as the reference configuration, B as the current conhguration and leaving B as 
the stress-free configuration, see Figure the initial stress becomes cr, the Cauchy stress 
becomes r, and Eq. (24) becomes 


T = p,C-^ + F-^aF-^ - pri, 


(25) 


where p^ is given by replacing r for a in Eq. (18)i, and Pr is an undetermined scalar 


Substituting (24) in (25) we find the connection 

{Pt-p)C~^ = {Pr-°P)I- 
As this equation must hold for every C we conclude that 

Pr = P and Pr = p- 


(26) 


(27) 


Note that, as is expected of a neo-Hookean material, the above equations do not determine 
p in terms of r or 

3.1 Plane Strain 


The free energy density (23) is simplified when the residual strain has only planar com¬ 


ponents. Accordingly, let us assume Bis = -B 23 = 0 and Bss = 1, which substituted in 


Eq. (16) results in ris = r 23 = 0 and tss = p — p. We let Bp, Cp, Bp, Ip, ap and rp 


be B, C, B, I, a and r restricted to the (xi,X 2 ) plane respectively. From equation (16) 
we get 

pBp = Tp + pIp. (28) 


To obtain p in terms of rp we take the determinant of each side of Eq. (28) giving 
det (/iSp ] = det(Tp-|-p/p) =i> = p -t-ptrrp-|-detrp 


(29) 


where we used det{pB) = p^. We solve the above for p to get 

P± = ± ^v^-4detr -h (trr)^ Ap^, 


(30) 


where —4detr -|- (trr)^ -|- Ap^ = (ri — T 2 Y + terms of the eigenvalues of r. The 

solution p = p_ should be discarded due to the following: p_ = —p when rp = 0 , and 


Eq. (28) shows that p should be positive when rp = 0, so for rp = 0 the only viable 
solution is p = p^. Further, as we expect p to be continuous in rp for every rp G 
and p_ 7 ^ p_,_ for every rp, we should discard the solution p_. 

For simplicity we assume Cs 2 = Csi = 0 and then Eq. (15) reads 


T = I (tr(BpCp) - 2 j + |(C33 - 1), 


(31) 


and substituting Bp from Eq. (28) we get 


^ = ^tr(rpC'p) -F ^ptrCp - p + ^{Css - 1 ), 


(32) 
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where p = p_^ given by Eq. (30). The Cauchy stress Eq. 0 becomes, 


ap = pBp - pIp + FpTpFp, 


( 33 ) 


and £733 = p — p + T 53 = p — p, where we have used that T 33 = p — p. Eor F = I we have 
<j = T and p = p. 


4 The Initial Stress Symmetry Equations 


as a set of scalar 


In this Section we explained how to express the ISS condition 
equations that relates the free energy density T; p and Pr from Eq. ( 1^ ; and the invariants 
of T and C. Let us first consider a simple example, namely assuming T = IiJi/2. One 
can check that this strain energy satisfies the compatibilty equations Q with p = trr, 
but does it satisfy the ISS (|12[)? We can use the stress in the form 0 to write 


<T = ^{F, T,p) = tr(TC)S — pi + tT{C)FTF'^, 
and then from ISS we have that 

r = ^{F-^,(T,pr) = - Pri + tr {B-^) 


(34) 


(35) 


To check that both Eqs. (34) and (35) hold for every r and F we substitute cr into Eq. (35), 


which after some rearranging becomes 

(1 - hh)T = (3Jl - 2pl2 + hIrl)C-^ + {hJl - Pr)I, (36) 

where we have used that tr(F~^) = I 2 , which can be shown by applying the Cayley- 


Hamilton theorem to B with 13 = 1 (due to incompressibility). In order to satisfy Eq. (36) 
for every F and r, the coefficients of r, C~^ and I must all be identically zero (this is 
shown rigorously in Appendix [B|) . For the coefficient of r to be zero, the identity I 1 I 2 = 1 
must hold for every F, which is obviously impossible: so we conclude that T = IiJi/2 
does not correctly furnish the Cauchy stress for every reference configuration. 

In Appendix we reduce ISS to nine scalar equations, following a procedure similar 
to the one above. They are compactly written as 




J 2 


+ T’jS'P* I = 0. 


'P.',. (ojl'i'i'A + + Qg'l'Pi + QjS'P*) = 0. 

where n sums over 0,1, 2,3,4, T := 1, 

^4 = ^h{F,r), = ^j^{F,t), 

:= ^i,{F-\a), := ^UF-\a), 

for A; G {1, 2} and m G {1, 2, 3,4} with a given by Q, 


b = 


(37) 

(38) 


(39) 

(40) 


to 




P 

-2^ 


p/2 


I 2 P - 2^p - 2 / 1 'I '/2 

Pt 

, = 4 

-'I'U 

II 

oT 

24-/2 

0 

-2$j, 

0 

0 


2l2'I'j2 


-24-j, 

1 




-44-/2 


— 0 


'{ 2 } 


'{ 1 } 


'{ 2 } 


(41) 


(42) 
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while all over matrices are given in Appendix]^ Eqs (37) and (38) represent 6 and 3 scalar 
equations, respectively, with the unknowns being the functions 'k, p and Pr- 

We now look at special cases of \k which simplify the ISS equations and lead to more 
practical representations for free energy density. 


4.1 Free energy independent of J 3 and J 4 


When 4' is independent of J 3 and J 4 , Eq. (38) is identically zero and only the first three 


terms of Eq. (37) are non-zero. The fourth scalar equation of (37) becomes —= 0 
which is true, for every F and r, only when T does not depend on J 2 or does not depend 
on Ji. We will investigate the first case in more detail below. 


Assuming T independent of J 2 


If T does not depend on J 2 , J 3 and J 4 then Eq. (37) reduces to 






h 


+ 






= 0, =4'q, 


(43) 


where we have used Eq. (43)i to derive the other three equations. We can check if these 
equations are consistent with the compatibility Eqs. ([^ by letting F = I, a = t and 
Pr = p = p- This results in = ±1/2, 2'kp = ±p and 'k/ 2 (l ± 1) = 0, which only 
satisfies Eqs. ([^ if Tjj = 1/2, 2Tp = p and T/j = 0 for F = /. 

Let us apply Eqs. (43) to the following example 


4 /= 2MA-3)±2(Ji-/n), 


(44) 


used by Merodio et al.(2013). In this case, Eqs. (43)i and ( 43)2 are satisfied while 


Eqs. ( 43)3 and ( 43)4 become p = Pr = p- However, restricting p to be constant, which is 


normally determined through boundary conditions, would result in physically unexpected 
behaviour. This is best seen with an example. 


Example 1 For the free energy density expressed in Eq. (44) the Cauchy stress from 
Eg.(§ becomes 

a = pF^F -Ip + F^tF. (45) 

Let us consider a cube in the reference configuration, subject to the residual stress Tij = 
Tidij, aligned with the axes of the cube. For the current configuration we impose two 
clamped conditions that fix Fn = Ai and F 33 = A 3 , then F 22 = A 2 = (AiA 3 )“^ due 
to incompressibility, with all the other components of F being zero. With this imposed 
deformation the Cauchy stress is diagonal with components aij = dijai. We can also 
prescribe the stress 02 , as this will not alter Ai and A 3 because they are kept fixed. We 
choose £72 = T 2 , to ensure compatibility with the reference configuration when Ai = A 3 = 1, 
which results in 

T2 = pXl -p + XIt2 . 

Clearly it is not possible for p to satisfy the above and be fixed at p = p so that ISS is 


satisfied. This means that the free energy density Eq. (44) either results in non-physically 


behaviour, or it does not satisfy ISS; which implies that the constitutive relation (45) does 
not hold for every reference configuration. 
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The choice of T should satisfy Eqs. (43) without restricting the deformation F, initial 
stress T or j?, otherwise the material will likely exhibit non-physical behaviour. Below we 
deduce a a free energy density T for a compressible model, inspired by the neo-Hookean 
example in Section that satisfies ISS without any unphysical restrictions. 

For a compressible material we substitute p with and Pr with —2/^^'!'^ in the 

ISS Eqs. (43). We will find under what conditions does the following free energy density 




(46) 


satisfy ISS (43) and stress compatibility 0 : where / and g are arbitrary scalar functions. 
The Cauchy stress ([^ becomes 


with 


= g{T, h)I^ + 2/3 4 /,3/ + Jg 


= \fhi.r,h) + \gi,{T,h)hI-^/^ - [g{T,h)h + Ji] ■ 


(47) 


(48) 


2 '' ^ 2 '" o g 

The ISS Eq. ([43l)i is satisfied as 'kjj'I'j^ = ( 1 / 2)13 ^^^( 1 / 2 ) 73 ^^ = 1/4 and stress com¬ 
patibility 02 is satisfied as 24'= 1 when F = I. By substituting p for — 2 / 34 ' 7 - 3 , the 
remaining compatibility Eq. ([^1 becomes 


g{T,I^) = —24 '73 evaluated at F = I. 


(49) 


Eq. (49) is satisfied if we set 1) = trT/3 — 3gi^{T, 1), then Eq. (49) and Eq. ( [47| ) 

together imply that cr = r when F = I. 


The ISS Eqs. ( 43)3 and ( 43)4 now read respectively 

-2I^^^l=g{T,h), 


- 2l34'73 = 


3 -') 


and 


(50) 


for every F and r, where fi^ and gj^ are respectively the partial derivatives of / and g 
with respect to I3. Note that the above two equations are equivalent as Eq. (50)i becomes 
Eq. (50)2 when we substitute F for F~^ and r for cr. For this reason it is sufficient to 
satisfy Eq. (50)i. If E' = I, Eq. (50)i is satisfied as it reduces to Eq. (49). The rest of this 
section will focus on showing that Eq. (50)i is satisfied for every F and r. 

By rearranging Eq. (47) and taking the determinant on either side we get 

det (cr - 2 / 34 ' 73 /) = det (i~^^^F^tF + g{T, b) = detEdet(r + g{T,h)I), 


(51) 


and by using incompressibility det S = 1 we find that 

det(cr - 2Is^j^I) = det(T g{T, 

If we choose g{T,I-i) such that 

det(r -h g{T, h)I) = k, 


(52) 


where fc is a constant, then the right handside of Eq. (51) is also equal to k, thus cr must 
satisfy 


/3/^det(cr-2/34'73/) = k. 


By swapping r for cr and F for F ^ Eq. ( |52[ ) becomes 

Ig'^^ det(cr -|- g{cr, = k. 


(53) 


(54) 
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For g{T,I‘i) to satisfy Eq. (52) there are three possible solutions given by replacing g? 
with and g{T,Is) with p in ( 18 )i,( 18)2 and ([ 1 ^ 3 , which we respectively denote by 

gi{T, I3), 5'2( t , I3) and gsir, I3). We highlighted in Sectionj^that only g{T, I 3 ) = 51 (r, I 3 ) 
is a physically viable solution. So if we choose g{T,l 3 ) = gi{T,l 3 ), then for Eq. (50)i to 
be satisfied we need 


( 7 i(cr, /g ^) = — 2/3'I'/3 for every F and r. 


(55) 


We will show that the above is a consequence of Eqs. (53), (54) and stress compatibil- 


conclude from Eqs. (53) and (54) that — 2/3 lb/g will equal either gi{ 


^), 52 ( 0 -,/g or 


ity (49). Due to Eq. (49), we know that initially for F = I Eq. (55) is true. We can also 

( 73 (<T,/g'^) for every F. Since — 2 / 3 'b /3 should be a continuous function of F and r, it can 
not change from gi{cr,I^^^ ^ r-i> 

T-i\ -L „ r-i 


) to another solution gk{a,I^ ) becanse it can be shown that 
gk{(T,I^^) / gi{a,I^^) for every cr and I 3 . Eor this reason if g'i(r, 1) = — 2 'b /3 for F = I, 
then gi{cr,I^^) = — 2 / 3 '!' 7-3 for every F. 


In conclnsion, the model (46) satisfies ISS as long as g{T, I3) = gi{T, I3), where 5 i ( t , I3) 
is given by replacing with in (18)i. 


5 Homogeneous stress gradient in a hollow cylinder 


It is now well acknowledged that residual stresses in living materials are vital to maintain 
ideal mechanical conditions. When an external load changes, growth and remodelling 
will alter the residual stress to best adapt to the new load. An excellent example is how 
arteries remodel in response to the internal pressnre. The residual stress in arteries is 



that if a given tissue grows in response to stress, then homoeostasis is only possible if the 
tissue is under similar stress conditions throughout. So by minimizing the stress gradient 
we are selecting the most homogeneous stress possible. We will also consider a simplified 
artery with only one layer and no shear stress applied to the interior of the artery. 

Our workflow is to first choose an optimal Cauchy stress field cr, and then the to 
derive the residual stress from the ISS equation r = ^{F~^,a,pr). So first we calculate 
the Cauchy stress with a near homogeneous circumferential and radial components in 
Section |5.2[ then we find the residual stress that supports this optimal Cauchy stress in 


Section 5.3 In Section 5.2 we also show that the optimal Cauchy stress has a simple 


asymptotic formula. 

We finally compare the results with the corresponding ones obtained using the opening 


angle method (Chuong and Eung, 1986). Pasquale: cute the rest. 


5.1 Plain Strain Cylinder 

To describe the arterial wall we use cylindrical coordinates (r, 6 , z) with the components 
of a written in terms of unit basis vectors. We assume that the arterial wall retains its 
cylindrical symmetry when the internal pressure is removed and that there is no shear 
stress at the inner wall. These simplifications imply that when a pressure is applied in the 
cylinder the resulting deformation is an inflation. 


12 

































Due to radial symmetry, the Cauchy stress cr is independent of 9 and z, so in the 
absence of body forces the equilibrium equations reduce to 


d / 2 
d 


d 

= 0 , —{razr) = 0, 


— {rarr) - (Tee = 0 , 

with the boundary conditions 

= 0 , azr = 0 and arr = —P for r = a, 
O'er = 0 , azr = 0 and arr = 0 for r = 6 , 


(56) 

(57) 


(58) 

(59) 


where a and b are the inner and outer radius of the loaded artery, respectively. The 


equilibrium Eqs. (56) together with the boundary conditions lead to aor = 0 and azr = 0 
for all r. As neither aze nor azz appear in the equations of equilibrium, they are only 
restricted by the constitutive choice and boundary conditions on the cross-section of the 
artery. We use this degree of freedom to assume there is no axial torsion aze = 0. Note 
that for the constitutive choice Eq. Tzz can be chosen so that azz is constant, while 
for the simpler choice of plain strain Eq.(33) (with azz = < 733 ) azz is determined by p. 


Finally, the remaining equation of equilibrium we need to enforce is Eq. (57). 


5.2 Minimal stress gradient 

Here we develop a method to minimize the stress gradient fields a'rr and a^g, where the 
prime denotes differentiation in r. We make no assumptions about any reference configu¬ 
ration nor do we make any constitutive choice. We only make use of the assumptions from 
the section above which result in both arr and aee being independent of the coordinates 
6 and z; aor = azr = Oze = 0 ; and for simplicity we will not consider azz- 

Once aee and arr are determined, we use ISS Eq.([IT|) to write the residual stresses r©© 


and trr in terms of aee, Orr and the deformation gradient in Section 5.3 


Using the equilibrium Eq. (57) we write aee in terms of arr as 

j 


oee — Orr ~\~ far 


age = ‘la'rr + ra: 


(60) 


As there are aortas and veins of many different sizes, for the sake of generality let us 
introduce the following dimensionless variables: 


Q = 


which result in 


and 


^ rr ^ 


r — a 
b — a 

P 


and c,{p) = 


o'ee — 


b — a' 
P 


b — a 

Our aim is to minimize the stress gradient density 


Orr(r) 

p ’ 
p 

{b — a)2 ’ 

(2? +ip + ay) ■ 


— r-" 

^ rr ^ 


(61) 

(62) 

(63) 


1 


rb r 


b — a 


[cr, 


rr) + i^ee)" 


p2 

dr = -^o? 




(?')^-h (2?'-h (p-Fa)?")^ dg (64) 


with the constraint 


?( 1 ) = 0 and <^( 0 ) = —1 


[ q'dg = 1, 

Jo 


(65) 
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where we have introduced the quantity a = a(b — a) 


Using the calculus of variations (Gregory and Lin, 1992) we find that must satisfy 
the Euler-Lagrange equation 


where 


df d df cz in u A n f n 1 

- deW’= ^^ ^ = »tai> = o,i, 

/ = (?0 ■*" + (^^ + «) 1 


( 66 ) 


(67) 


and A is a Lagrange multiplier due to the constraint ( |65| ). The solution to (66) is given 
by 


. » LS-lI I 1 

, A A , r— .a 2 ■^ 2 — (1+q;) 2 ^2 ^v^ 4 _l 

? = A - e(^- 3)--4t^(« + e)~ 


6 6 


+ 


2A 


Q!)^bi3 

a 2 A — (1 + a)^ 


n/ 13 1 . , yi3 1 1 

' -L 


3(vT 3 — 3) a'Tbl _ (1 -)-q/^a/Is y a + ^ 


We determine A from the constraint (65) to be 


%/13 I 1 

a{l + Q^)^ 2+2 


( 68 ) 


A 

6 


x/l3 


1 - (\/l3- 3) 


a 


us I 1 , via I 1 

2+2— (1+ a) 2+2 vU3 

(a + £>) 2 - 


a'Tbi — _)_ q;)\/I3 

v/13 


+ 


a 2 


1 /- ' \/l3 1 1 

2 — (1 + a) 2 '-2 


VTs — 3 a'/bi _ + q;)’7I3 


s AH- 

a(l + a) ^ 2 ■ 

a + Q 


dg) \ (69) 


Realistic values for a = a{b — a)~^ can be obtained from in vivo measurements of the 
lumen radius (inner radius) a and total artery thickness b — a. For the descending thoracic 
aorta of healthy men around 51 years old the mean lumen radius is 20 mm (Stefanadis 


et ah, 1995) and the mean thickness is 1.4mm (Mensel et ah, 2014). These estimates give 


a ^ ~ 0.07 making it appropriate to expand and A as a power series in a which gives 


^ (l + (1 - 2p)a-^ + ^ (9p2 - 6p - 1) + 0{a-^) , (70) 

A = 2 + —ot ^ + 0(a , (71) 

o 

or simply 

?' = 1 + (1 - 2p)a-^ + p(3p - 2)a-2 + ©(a”^) . (72) 

We can substitute q' in Eqs. ( [60| ), (61) and (63) to obtain 

Urr = -P{^ - /o)(l - + P^oT'^) + 0{a~^) , (73) 

aee = Pa{l + a~^p‘^{Ap - 3)) + 0{a~^) . (74) 


We plot the above Cauchy stress for a~^ from 0.05 to 0.1, which is within a physiological 
range, in Figure]^ A graph for arr/P would show all the curves bunched along the same 
straight line, so instead we have depicted the curves for arrjP — P- Figure [^reveals that 
as the relative thickness of the arterial wall increases (shading from blue towards red), the 
circumferential stress aoe decreases while the radial stress cr^r both increases and becomes 
less homogeneous. 
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Figure 2: The plots of the dimensionless Cauchy stress against the dimensionaless radius 
p, where the color shades from blue to red as a~^ (the wall thickness divided by the inner 
radius) goes from 0.05 to 0.1. 


A well established method to quantify the residual stresses within arteries is the open¬ 


ing angle method (Chuong and Fung, 1986). For a neo-Hookean material (13) we will use 


the opening angle method to minimize the circumferential stress component 

1 /daee 

b - a Ja \ dr ) 


(75) 


in terms of the opening angle (p, restricted to the boundary conditions (58) and (59). The 


only two parameters with a unit of time, and mass, are P and p, so by rewriting P = P^p 
the stress will not depend on p. The results for Pq = 0.05, Pq = 0.2 and Pq = 0.35 are 


compared with the optimal stress Eqs. (73) and (74) with parameters for the descending 


thoracic aorta in Figure]^ We can see that aee for the opening angle method converges 
to the optimal stress aee as P/p tends to zero, while the plots for cr^r all overlap. Note 
however that azz for the opening angle method is not necessarily homogeneous in r. 


PeelP aJP-p 



Figure 3: The graphs of the dimensionless Cauchy stress against the dimensionless ra¬ 
dius p for the opening angle method (dashed curves) and the optimal stress Eqs. (73) 
and ([^ (red solid curves), for which we used parameters for the descending thoracic 
aorta: a = 20 mm, b = 21.4 mm, a~^ = 0.07 and = 1. As the dashed curves shade from 
yellow to green P/p goes through the values 0.05, 0.2 and 0.35 with optimal opening angle 
(j) = 65.9°, 204.2° and 261.723°; stress free reference inner radius 20.3mm, 25.2 mm and 
31.1mm; and stress free reference outer radius 22 mm, 27.8 mm and 34.4 mm respectively. 
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We can now invoke Finite Elasticity to derive the residual stress r in the unloaded 
state necessary to sustain the optimal homogeneous Cauchy stress, and then compare our 
results with the residual stress predicted by the opening angle method. To do so, we let 
T be a function of both stress and strain and use the initial stress symmetry (ISS), in the 


form of Eq. (11), to determine r as a function of a and F. 


5.3 Residual stress 

Here we use Einite Elasticity to connect the current state with the unloaded state. Since 
the cylindrical symmetry of the artery is maintained when the internal pressure is re¬ 
moved, we describe the unloaded state with the cylindrical coordinates (i?, 0 ,Z), then 
the deformation gradient for the unit basis vectors e^, eg, e^, Er, Eq and Ez becomes 


dr 


Er + ® Eq + Xz^z ® Ez-, 


(76) 


where we let be a constant. Assuming the material is incompressible we get 

dr r . 

1 ^ 

Let the reference configuration be a hollow cylinder with inner and outer radius A and H, 
respectively, so that 


(77) 


r = ^(IP - ^ ^ = 


-1 
2 : 5 


(78) 


which we will use to replace dr/dR wherever it appears. We will assume that the residual 


stress T is homogeneous in 0 and Z, as the deformation,the Cauchy stress from Eqs. (73) 


and (74), the boundary conditions (80) and (81) are all homogeneous in 0 and Z. 


The equilibrium equations now reduce to 

tqq = {Rtrr)' iox A< R< B, 

with the zero-traction boundary conditions 

trr = 0 for R = A, 
trr = 0 for R = B. 


(79) 

(80) 
(81) 


After making a constitutive choice for 'L, the residual stresses trr and tqq will be com¬ 
pletely determined from F, a and Pr due to ISS 0 

T independent of I2, J21 J3 and J4 

To simplify we assume T independent of I 2 , J 2 , J 3 and J 4 , so the Cauchy stress (|^ 
becomes 

a = 2 ^rFF'^ - pi + 2 ^,rFtF'^, (82) 

and, from ISS (12), we can swap r, F and p respectively with a, F~^ and Pr to get 

T = 2'ff^j^F-^F-'^ - pri + 2T5^F^icrF“^. (83) 


Substituting F from Eq. (76) into the above equation we arrive at 

2 

trr = 2 A 2 ^ ( TJ ^ + ^\arr) - Pr, 
r2 

'^ 0 © = 2 ^( 4 '?! + ^RO-0e) - Pr- 


(84) 

(85) 
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We remark that tzz has not been derived since we have not taken into consideration Ozz 
in Section 

To compare unloaded arteries of different sizes we write r and R in terms of the 
dimensionless radius p, 


r = {b — a)p + a and R= \/ + \z[{a + {h — a)pY — o?\, 


which we use to rewrite the ODE (79) in the form 


jjdTRR jjdTRR dp 

R + trr — Tee = R —+ trr — tqq 


dR 

dTRR 

dp 


R 


dR 

dp 


-1 


dp dR 

+ Trr - Tee = 0, 


which implies that 


dTRR I dR 


( 86 ) 

(87) 

( 88 ) 

(89) 


Integrating both sides in p and using the unloaded boundary condition (80) and (81) we 
reach 

f^TQQ-TRRdR f^Tee-TRRdR 

trr= / - 5 - —dp = 0 and / - - - —dp = 0. (90) 

Jo R dp Jq R dp 


Eq. (90)2 is independent of Pr, and can be used to solve for one of the parameters P, a, 


b, A, Xz or p, note that B depends on the other parameters because Eq. (78) evaluated at 
R = B gives B"^ = -\- Xz{b‘^ — a^). The only two parameters with a unit of time or mass 


are P and p, so by rewriting P = Pop the Eqs.(90) will not depend on p. To illustrate, we 


adopt the neo-Hookean material model (32), which results in = 1/2 and 2^^^ = p{cr) 


where p is given by substituting r for cr restricted to (r, 9) in p_^ of Eq. (30). 

To calculate the residual stress that supports the optimal Cauchy stress, we substitute 
aee and arr Eqs. (73) and (74) into re© and trr Eqs. (84) and (85). We can then use 
Eq. (90)2 to determine the unloaded geometry from the loaded geometry. To illustrate 


we take the parameters for the descending thoracic aorta, as used in the previous section, 
a = 20 mm, b = 1.4 mm, and = 1. We can then determine the unloaded inner radius A 
for different values of Pq, which is shown in Figure Surprisingly the unloaded geometry 
given by the opening angle method is approximately the same as shown in Figure 

We can also investigate the residual stress in the unloaded state. Once A is determined 
from Eq. (90)2 we can use Eq. (90)i to determine trr, and then r©© from Eq. 


To 

compare with the opening angle method we choose P/p = 0.05, 0.2 and 0.35, which all 
give a reasonable unloaded geometry, see Figure]^ The results are shown in Figure 

Biological tissues are very energy efficient. It is a common line of reasoning that 
biological systems adapt so as to minimize their potential energy. So it is possible that 
biological materials remodel their residual stress to lower their potential energy. Figure 
below shows the free energy density T for the opening angle method and for the model 


with the optimal stresses Eqs. (73) and (74). In all cases T is approximately a straight line 
that descreases as p increases away from the pressured boundary p = 0 (r = a). Figure]^ 
shows the difference between T/P from the opening angle method minus T/P from the 
model (32), which we denote as AT/P. The integral of AT/P over p G [0, 1] is positive 


if the opening angle method has on average a larger free energy density than the model 
(32), with the optimal stresses Eqs. (73) and (74). We found that as P/p increased, so 


did the total free energy in the cylinder cross-section AT/P(6 — a){p{b — a) + a)dp; for 
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Figure 4: The loaded artery wall is illustrated by the red curves, while the unloaded artery 
wall is illustrated by the blue curves. The y-axis illustrates the position of the loaded and 
unloaded walls measured from the center of the artery. The stress in the loaded geometry 
is taken to be the optimal Cauchy stress (73) and (74), and we assume a residually stressed 
neo-Hookean material (32). 


Tqq/P ^rr/P 



Figure 5: shows the dimensionless residual stress against the dimensionless radius p. The 


dashed curves are from the opening angle method and the solid curves are from Eq. (79) 


together with the optimal stresses Eqs. (73) and (74). As P/p. goes through the values 


0.05, 0.2 and 0.35, the dashed curves shade from yellow to green and the solid curves shade 
from blue to red. The unloaded geometry for both methods is approximately the same 
and shown in Figure]^ The other parameters are a = 20mm, b = 21.4mm, a~^ = 0.07 
and Az = 1. 


instance, this integral evaluates to —1.6410 ^ for P/p = 0.05, 2.3710 ^ for P/p = 0.2 and 
4.2510-3 for P/p = 0.35. 

Essentially we can see from the Figures and that the optimal stresses Eqs. (73) 
and (74) with the neo-Hookean material (|32|) produce a more homogeneous stress and lower 


free energy than the opening angle method, though the two methods gave very similar 
results. The advantage in using the ISS (12) for a homoeostasis hypothesis for the Cauchy 
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stress cr is two fold. First is the freedom to choose cr so as to satisfy the homoeostasis 
hypothesis, which with ISS becomes a separate step from choosing a constitutive equation. 
Second is that it can be relatively straightforward to use ISS (12) to quantify the residual 
stress needed to support the chosen cr. 





Figure 6: (a) compares the dimensionless strain energy density ^/P versus the dimension¬ 
less radius p resulting from the opening angle method, yellow (green) dashed curves, with 
the model (32) whose residual stress maintains the optimal stress Eqs. (73) and (74), the 
red (blue) solid curves. The parameters used are given in Figure]^ (b) shows T/P from 
the opening angle method minus T/P from the model (|32[), which we denote as AT/P. 


6 Conclusions 


In order to quantify the initial stress within a solid using non-destructive experimental 
techniques, it is convenient to write the free energy density T = 'I'(P,r) in terms of the 
deformation gradient F and the initial stress r. 

In this article we presented a new constitutive condition, the initial stress symmetry, 
that aids in proposing suitable constitutive relations for T = 'I'(P,r) by providing nine 


scalar equations, see Eqs. (37) and (38). One immediate result is that guessing a functional 
dependence for 'I' = 'I'(P, r) is not a trivial task. In fact all choices for 'k(P,r) used in 
the literature do not satisfy ISS. Conversely, using ISS, we proposed two simple choices for 


\J/ = T(P,r), one incompressible Eq. (23) and one compressible Eq. (46), which include 
a generalization for an initially stressed neo-Hookean material. 

One consequence of ISS is that the initial stress can be derived from the Cauchy stress. 
Furthermore, ISS suggests that it is possible to first choose the Cauchy stress and, second, 
to make a constitutive choice and then use Eq. 0 to determine the corresponding initial 
stress. We illustrated this method in Sectionj^for a simplified arterial wall, where we chose 
the ideal Cauchy stress by using a minimal stress gradient hypothesis as the homeostatic 
condition, and then we calculated the optimal residual stresses of the unloaded remodelled 
artery. 

Since initial stresses are widespread in both inert and living matter, the proposed con¬ 
stitutive theory can be used in many applications towards a non-destructive quantification 
of initial tensions in solids. Euture research include determining of a wider class of material 
behaviours, e.g. including a natural material anisotropy, and linking the initial stress to 
elastic wave speeds. 
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Appendices 

A Stress-free configuration implies ISS 

Referring to Figure the virtual stress-free configuration guarantees that for any r and 

O 

F there exists a F such that 

a = 'd{B,p) and t = 'd{B,p). (A.l) 


^ _ 1 ~ _ T' ® 

where B = F~ BF~ ; p and p are a scalar fields where p depends on the boundary 

~ o o o 

conditions on B and B, while p depends on the boundary conditions on B and B] and 

54 ' 


■d{AA\q) := A 


d(A^A) 


A^ - ql. 


By setting r = 0 in Ec 


Johnson and Hogei (1995) 


Q we can see that the right side depends only on AA?" and q. 


demonstrated that it is possible to determine B from any 

O 

given T, together with appropriate boundary conditions on B, in order to reach useful 
constitutive equations such as 


-^{FBF^,p) = ^{F,t,p) for any F and S. 


(A.2) 
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As Eq. (A. 1)2 is valid for any F and B, we can substitute 


F for F 


B for B 


(A.3) 


and swap the boundary conditions on B and B (these substitutions effectively swap B 


for B), so that Eq. (A.1)2 becomes cr = 'd{B,p), where we have assumed that the above 


together with the boundary conditions on B has a as the unique solution. Analogously, 


Eq. (A.l)i becomes r = i9(B,p). 


This means that when we make the substitutions (A.3) in Eq. (A.2) we should also 
swap T and cr, so that 

= ^{F-^,a,p), (A.4) 


The left-hand side is simply r and the right-hand side is r given by ISS (12). Finally, the 


above condition is identically true, since this result is valid for every F and B. Moreover, 


we assumed that B can be determined from r, so Eq. (A.4) holds for every F and r. 


B Reducing ISS to nine scalar equations 

Deriving the scalar equations equivalent to ISS requires lengthy calculations. First we use 
the stress expansion and ISS ([T^ to write the residual stress as 


T = -pri + 2Tf^C"^ -F - C-^) + 2T5^F~VF-^ -F 2'^’}^F-^(7‘^F-^ 

+ 2'if%F-^{(7B-^ -F B-^ct)F-^ + 2'i'’}^F-^{a'^B-^ + B-^a‘^)F-^ (B.5) 

where I 2 = tr(C'~^) from Cayley-Hamilton with detC = 1. The above, together with (§ 
for cr, must hold for every r and C. To make use of this restriction, we substitute cr 
from ([^ into Eq. (B.5) to write an equation in terms of only t, C, p and p-r, compactly 
written as 


= 0 , 


(B.6) 


where we adopt the convention that repeated indices in a single term implies summation 
over all the values of the index, with i,k,n G {—2, —1, 0,1, 2} and j,m G {0,1,2}. The 
aijkS can be calculated with a computer algebra system and are too cumbersome to 
reproduce here, though it is important to remember that the atjkS are functions of p 
and Pt] the invariants (§, @ and (©; and the derivatives of T with respect to the 
invariants ([^ and ([^. Due to the Eqs. (B.5) and ^ being symmetric, we have the 
symmetries Otnjkmi ^nmkji- 

Our aim here is to use Eq. (B.6) to obtain the corresponding scalar equations that 


hold for any choice of T and the invariants ([^, (Q and (§. To acheive this, we need to 
understand how we can vary C and r while keeping the invariants hxed. We can then 


apply these variations to Eq. (B.6), which holds true for any C and r, to reach scalar 


equations in terms of the invariants and T. 

So first, to keep the invariants (© and @ fixed, we have to keep the eigenvalues of 
both C and r hxed. We can assume that the eigenvalues Af, A| and A| of C are all 
different, as this is true for most any deformation and so we do not loose generality. The 
same assumption can be made for n, T 2 and rs, the eigenvalues of r. 
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The mixed invariants involve the eigenvectors of C and r. So we write the eigen 
decompositions 

3 3 

C = '^^XfVi^Vi and t = '^^Tjtj i^tj, (B-7) 

i=l j=l 


where the Vj’s and tj’s are respectively the eigenvectors of C and r. Using these decom¬ 
positions we see that 


Ji = tr(TC) = . tjf, J 2 = tr(TC2) = TjXf{Vi ■ tjf, (B. 8 ) 

J 3 = tr(T2C) = r^Af (F, • tj)\ J 4 = = r|A^(F, • (B.9) 

where summation over i and j is implied. Setting the above four terms to be constant can 
be seen as four independent equations for the Fj’s and tj’s. Being unit eigenvectors of 
symmetric matrices, they must also satisfy the twelve equations 


F i ' F fYi — ^im and tj • — ^jn: 


(B.IO) 


for i,j, m,n = 1, 2 ,3. 

So in total, the F^’s and tj’s have to satisfy 16 equations for the 18 unknown compo¬ 
nents of Fi, F 2 , F 3 , ti, t 2 and < 3 . This leaves us with two degrees of freedom which we 
can carefully use to reduce Eq. (B. 6 ). For example, let Fi -ti = F 2 -ii = 0, which implies 
that ti = TFs and so ti is an eigenvector of both C and r. By taking the dot product 


of Eq. (B. 6 ) with ti on the left and right side we conclude that 




j-\-m^{i-\-k-\-n) _ 


= 0 . 


(B.ll) 


Now, in the same way, we could have used the degrees of freedom to choose V i ■ ±2 = 
F 2 • *2 = 0 and concluded that 


j-\-m \ 

(^ijkmn'^2 


j+m \2{i+k+n) _ 


= 0 . 


and similarily choosing Fi • ts = F 2 • is = 0 we conclude that 

= 0 . 


(B.12) 


(B.13) 


We can rewrite the three above equations by using the characteristic polynomial of A^ and 
of Tp, 


A3 — /1A3 — /2A3 -|- 1, '^3 — — It 2 T 3 + It 3 j 


(B.14) 


to replace powers of A 3 higher than 4 and lower than 0 with A^, A^ and the invariants of C, 


and the analogous for powers of ri, T 2 and In doing so Eqs. (B.ll), (B.12) and (B.13) 
respectively become 


/SijXfrl = 0 , PijXiW^ = 0 , PijXfr^^ = 0 , 




2 * J 


(B.15) 


for i,j = 0,1 and 2 , and where the /Ijj’s can be written in terms of aijkmn and the 
invariants of r and C. This is best done with a computer algebra system. 

As a reminder, the eigenvectors of C and r in each of the above three equations are 
most likely different. However, the eigenvectors only appear in the form of the mixed in¬ 
variants , which are the same for the three equations because the of Eqs. ( |B. 8 [ ) and ( |B.9[ ). 
Thus we may solve all the above three Eqs. (B.ll), (B.12) and (B.13) simultaneously. 
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The three equations state that ri, T 2 and T 3 are the roots of the second-order polynomial 
= 0. As these eigenvalues are assumed to be all different, this is only possible if 

A, Ai' = 0, for j = 0,1,2. (B.16) 


Finally, there was nothing special about V 3 and A 3 , using analogous arguments we can 
conclude that 

A,Af = 0, A,A2* = 0, A,Ai' = 0, fori = 0,1,2. (B.17) 

Again, Ai,A 2 and A 3 are assumed to be all different, so the only way that Ai, A 2 and A 3 
are the roots of a second-order polynomial is if the polynomials coefficients are all zero, 
leading to 

Pij = 0, for z,i = 0,1, 2. (B.18) 


This result suggest a simple way to reach expression for the /Ijj’s: we can substitute scalars 
A and r, in place of the tensors (B. 6 ), followed by repeatedly using the characteristic 
polynomials of both A and t as explained around Eq. (B.14), and then the /Ijj’s will be 
the coefficients of the multivariate polynomial in A and r. The results are explicity given 
by Eqs. (37) and (38) combined with Appendix 


C ISS matrices 


The matrices appearing in Eqs. (37) and (38) 


have the following expressions: 



2'1'f^ +2^-73 (p-/24'73) 


0 

p72 - 4^714'73 - 2/14-13 


0 

2/14-2^ + 2^23 - 24'7,(p - 2/2^73) 

= 8 

0 

0 

’ {1} ^ 

-24-71 

0 


2/24-7, + 24-73 

0 


j) - 2 /i4-7, - 2/24-73_ 


0 

0 

e){3} — 1 C 0 

{ 2 } “ p- - 2 / 2 '!'/, 

2(/i/2 - 1)^'/, + h(‘ih^h ~ P) 
JiP - 2(/i - + 2(1 “ 


Ir3{2^J, + 4/1 + Irl^Ji) 

4^3 

p{3} ^ A3(2 /i^J, + 4(/2 - /2)^+ hlrl^J,) 

{3} (T3 - IrlIr2)^Js ~ - 4/l/rt«'./, 

2I2It2^Ji + 4(/i/2 — 1)^241 ,/2 + (A 1 A 2 — lT3)h'^,h 
.4(^2 ~ Ii)It 2'^.h ~ 2 IiIt2'^Ji + (A 3 ~ AiA2)A4'j3 _ 


(C.19) 


(C.20) 


(C.21) 



32 


A3 (2 (A -/f) T72 - 
-A3 (2 /i^'j2 +^, 7 ,) 

As (2 ( 2 A /2 - A" - 1) + (A - A) ^Ji) 

A 2 (2 (A-A) 'E'j2+A4'j0 
A 2 (2 (A + A - AA) 4>j2 + (1 - AA)4>7,) 
A 2 (2 (A - 2 / 1/2 + 1) + (A “ A) ^Ji) 




-AiA3((A-A)^J4+A4'J3) 

-A1A3 (A4774 + 4 /,73) 

-A1A3 ((A - 2AA + 1) 4>j 4 + (A - A) 4773) 
(A1A2 - As) ((A “ A) 477, + A4I73) 

(A1A2 - As) ((a + A - AA) 4>74 + (1 - AA) 4 > 73 ) 
(A1A2 - As) ((A - 2/1/2 + 1 ) ^74 + (A - A) ^73) 


(C.22) 


24 




















p{4} _ 

^{0} 


= -4 


-4/ip'I'/ 2 - 4p'I'7j + + hp^ 

0 

0 

0 


p{4} _ 

P{1} - 16 


0 


0 

0 


0 

0 

24-/2 

p^^J' — 32 
, -f'|2} — 3/ 

0 

(M 

1 

p - 2/24-12 


2 (/24-I1 + 4-12) 

. -2^/3 . 


-2/i4-ii -2/24-I2 + p _ 


p{4} _ 

^{3} - 


-16 


4/^34'j, 

0 

/r3 (4/i4'j, -F/rf'I'j, -f 25'jJ 
-4/rt4'j2 

4/2/t24'J2 

.4'./3(/r3 ^ ItiIt2 ) ~ 4/i/t 24'J2 — 2/rt4'jj 


p{4} _ 

-^{4} - 


-64 


\ /rt (2/1^42 +^Jl) 1 


hlhi (-fl4-44 + 4 - 43 ) 

to 

to 


1 ^ 7 - 11 ^ 7-3 ^44 

hi (2 {II - h) 4-42 + d4-4i) 

-64 

/rl/r3 {{II - h) 4 I 44 + d4>43) 

“2/i/t24-42 - /t24-4i 

-{ItiIt 2 - hi) ih'^Ji + 4 - 43 ) 

/r2(2(/l/2 - 1)4-42 -//24-4 i) 


{hlh2 - hi) {{hh - 1)4-44 + ^ 2 ^ 43 ) 

h2(2 {h-lf) ^42-d4-4i). 


_{hlh2 - hi) {{h - If 4-44 - 


Q|i| = QI 2 I = ^{ 1 } = = Q{o} = n = 1, 2,3 and 4, 


*^{3} “ 


*^{4} - 


-1 

h 

-h 


— 

*^{3} “ 


0^2^ - 
<^{4} - 


-16 


<V{i} - 


n^^j' — 

‘^{3} “ 


-'4'Ji 

h^Ji 


V{2} - 


-32 


+ (/2 - I 2 ) 

(1 - hh)-^j, + {-hll +h + If) ^J2 

L {II - h) + {If - 2 / 2/1 + 1) ^J2 J 


- 2^jJrl - '^h^jjrl - 24 '/! -/ It2^J, 


+ 2/24'ji/tI + 4/i/24'j2-fTi - 44'j2-frt + 2h^i^ + ~ I2 It 2^.h 

lT2^j3h +P- 2/24'/2 + 4/ri/2'I'./2 - I^ih (44'./2-fi + 24'/i + 24'-//ri4'j3). 


(C.23) 

(C.24) 

(C.25) 

(C.26) 

(C.27) 

(C.28) 

(C.29) 

(C.30) 




r h-ii ] 


-h 

32 {If - hi) 4-44 

hll -h- II 

2/2/1 -1 - /? 

-/ 32 {If - hi) 4/j3 

hh - 1 
h - II 


16 


p - 2/24'/2 - 2/1 (4'/j-F/ri4'ji) + 4 (/2 -/?)/ri4'j2 
2 (/i/ 2 — 1)4//i -“ I 2 P ~ 2I.ri^ji + 2/2 (/24'/2 -I- /i/ti4'ji) — 4 {—hll 4- /i + if) /ti4'j2 , 
.2 {h - If + hp + 2(1 - /i/2)^/2 + 2 (/2 - If hi^j, - 4 (/? - 2 / 2/1 + 1) /ri4'j2. 


(C.31) 


p){4} _ 

^|1} — “16 

^{4} _ 

C^{3} — 16 


■ 0 ■ 

0 , 



= 64 


-4-41 - /i4'j2 
h'^Ji + {hh - 1)4/42 ) 

Xh-if ■^.h-h^.h. 


2 (4-12 - 2/,i4'42) 
p - 2/24'i2 + 4/2/^!4'j2 

4 - 43/21 - 24-41/,! - 4/i4'42/ri - 24'ii -/ /rf4'43 


(C.32) 

(C.33) 




■-1' 

1 

-h 

64 {If - hi) 4-43 

h 

1 

+ 

hh - 1 


-h 

1 

h - II 


32 


2 (/24'/i -I- 4/42 + hhl^Ji + 2(/i/2 - 1)/ti4'j2) 
-2 (4-41 -+-/t 1 (4-4i + 2 / 14 - 12 )) 

_p - 2/24-12 - 2/1 ( 4-/1 -/ /rl4-4i) + 4 {h - If) /^i4-j2 


(C.34) 
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